Thermodynamic properties of a simple model 
of like-charged attracting rods 
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We study the thermodynamic properties of a simple model for the possible mechanism of at- 
traction between like charged rodlike polyions inside a polyelectrolyte solution. We consider two 
polyions in parallel planes, with Z charges each, in a solution containing multivalent counterions 
of valence a. The model is solved exactly for Z < 13 for a general angle 9 between the rods and 
supposing that n counterions are condensed on each polyion. The free energy has two minima, one 
at 6 — (parallel rods) and another at 9 — tv /2 (perpendicular rods). In general, in situations where 
an attractive force develops at small distances between the centers of the polyions, the perpendicu- 
lar configuration has the lowest free energy at large distances, while at small distances the parallel 
configuration minimizes the free energy of the model. However, at low temperatures, a reentrant 
behavior is observed, such that the perpendicular configuration is the global minimum for both 
large and small distances, while the parallel configuration minimizes the free energy at intermediate 
distances. 
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I. INTRODUCTION 

It has been observed that like-charged macromolecules 
can attract each other in solutions containing multivalent 
counterions. This attraction manifests itself in in vitro 
formation of toroidal bundles of concentrated DNA 2|, 
similar to the one found in bacteriophage heads 
and the appearance of rodlike bundles of f-actin and to- 
bacco mosaic virus in the presence of multivalent coun- 
terions A number of models have been suggested 
to explain this curious phenomenon. The fundamental 
ingredient in all of these models is the role played by 
the condensed counterions (|J^. Thus, the attraction 
has been attributed to the correlations between the con- 
densed counterions on the two polyions . The math- 
ematical problem of how these correlations can be taken 
into account is highly non-trivial. Two approaches have 
been proposed. One relies on field theoretic methodol- 
ogy and uses what can be classified as a high tempera- 
ture expansion to account for the correlations between 
the condensed counterions ||^, while the second is a zero 
temperature approximation, in which the counterions are 
thought to form correlated Wigner crystals on the sur- 
faces of the two polyions |p|-pl| . Neither one of the ap- 
proaches is exact, although the zero temperature Wigner 
crystal approximation seems to be better at capturing 
the true nature of correlations . 

In a previous work we have introduced a simple model 
which has allowed us to study exactly the force between 
two parallel charged rods with a layer of condensed coun- 
terions M . In this paper we shall extend this work to al- 



low for a relative inclination between the two rods. This 
problem is of particular interest in kinetics of the bundle 
formation. Indeed, it has been observed that the bundles 
of stiff polyions have a characteristic size. This is quite 
surprising since the correlation induced attraction should 
favor the formation of infinitely thick bundles, after all 
the ionic crystals can grow to macroscopic sizes. A possi- 
ble explanation for the bundles not growing beyond some 
specific size can be found in the kinematics of bundle for- 
mation. The condensed counterions do not fully neutral- 
ize the charge of a polyion. For a rodlike polyion formed 
hy Z —q charges separated by distance b, placed in a sol- 
vent with dielectric constant D, in the presence of coun- 
terions of charge aq, the Manning criterion states 
that the number of condensed counterions on a polyion 
is n = (1 - l/a^)Z/a where ^ = q^/DbkBT. Thus, for 
the case of DNA with divalent counterions 88% of the 
DNA's charge is neutralized. It is easy to convince one- 
self that if the interaction between two rodlike molecules 
is repulsive, there will be a greater probability that they 
will be found perpendicular to one another. The corre- 
lation induced attraction between the DNA molecules is 
short ranged, what means that the electrostatic repulsion 
is dominant on large scales. Thus, the two polyions will 
in general repel one another. It is only when the two 
macromolecules come in a close contact that the thermal 
fluctuations might be able to overcome the free energy 
barrier between the perpendicular and the parallel con- 
figurations, allowing the correlation induced attraction to 
take over and the polyions to "bundle up" . It has been 
argued that the size of this free energy barrier scales with 
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the size of the bundle aheady formed [|^. Thus, there 
comes a point when the thermal fluctuations will not be 
able to overcome the free energy barrier. Motivated by 
this discussion we shall now proceed with the study of 
interactions between rotating like-charged rods. 



II. DEFINITION OF THE MODEL AND ITS 
SOLUTION 



We consider two rodlikc polyions, each one with Z 
charges —q with separation h between neighbor charges. 
On each polyion, n a-valent counterions, with charge aq, 
are condensed. The sites where the condensed counteri- 
ons are located are described by occupation variables oij , 



such that cr," 



1 if a counterion is condensed on the I'th 



monomer (z = 1, . . . , Z) of the j'th polyion (j = 1, 2) and 
<Jij = otherwise. When a counterion is condensed on 
a monomer, we assume that the only effect is the renor- 
malization of the local charge, from —q to —q -\- aq. The 
rods are located on two parallel planes separated by a 
distance d and the line joining the centers of the rods is 
supposed to be perpendicular to these planes. The angle 
between the directions of the rods is equal to 9, so that 
9 = corresponds to the case of parallel rods previously 
considered The definitions above are illustrated in 
Fig. |]. Notice that, for simplicity, we will consider only 
odd values oi Z, so that for nonzero values of 9 the dis- 
tance between the central charges of each rod vanishes 
as c? ^ 0. The polyions are placed in a uniform solvent 
whose dielectric constant is equal to D. For a given con- 
figuration {cr} of condensed counterions the Hamiltonian 
for the pair of polyions may then be written as 
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where the denominator is the distance between the sites 
and (i ,j ), the sum is restricted to 7^ (« ,j ), 

and 



x\i-i\, if j = j 



with X ~ d/b and 



x^+f! + B~2f,f.cos9, if j^j ; 



Z + l-2i 



(2) 



(3) 



The values of occupation variables obey the constraints 

Ez \ — ^ z 



FIG. 1. Two rodlike polyions with z — 7 charges each. 
The partition function of the model is given by 

g = ^' exp(-/3H) = cxp(-^i?), (4) 



where the prime denotes the constraint of fixed numbers 
of condensed counterions on each polyion. The adimcn- 
sional reduced Hamiltonian H is given by 



H = 



2 

,3') 



(1 - acry)(l - acr,,'y) 



(5) 



and ^ — f3q^ / Db is the Manning parameter [^. For a 
given counterion configuration {cr}, if a transformation 
a^j = 1 — cr.y is performed, the Hamiltonian changes as 

H{Z, a, n, X, 9, {cr}) = (a — 1)^H{Z, a ,n ,x,9,{a }), 

(6) 

where a — a/{a — 1) and n = Z — n. This relation 
leads to the following invariance property of the parti- 
tion function of the model 

Q{Z, n, ^, a, X, 9) — Q{Z, Z — n,[a — a/[a — l],x, 9), 

(7) 

and therefore we may restrict ourselves to n < Z/2 in 
the calculations. The free energy of the model is = 
— fc^rinQ, and the force between the rods is F = ~§d' 
It is then useful to define an adimcnsional force as 



Db^F 1 dQ 



(8) 



dx ' 

It should be stressed that the definition of / we use here 
is different from the one used in , since the earlier def- 
inition diverges at vanishing temperature. The sign was 
chosen in such a way that repulsive forces are positive. 
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To solve the model exactly, one may define activity 
variables 



exp 



exp 



i = l,...,Z-l 

ih j ~ ■ ■ ■ T 



(9) 
(10) 



the perpendicular {9 = tt/2) configuration of the rods. 
An example of this is shown in Fig. y, and one may note 
that, in general, the parallel configuration is stable for 
small distance x and the perpendicular one becomes sta- 
ble as X is increased. So, in what follows we will concen- 
trate our attention on the parallel and perpendicular rod 
configurations only. 



where dij — \Jx'^ + ff + ff, — 2/^/^/ cos 6. It may then 
be noticed that the partition function may be written as 



(11) 



where 



Nr. 



Z\ 



n\{Z -n)\ 



is the number of condensed countcrions configurations 
and the statistical weight of the i'th configuration is given 

by 



T[yr n 

j=l k,l=l 



(12) 



Uij and Viki being quadratic polynomials in a with integer 
coefficients. It is possible to generate these sets of inte- 
ger numbers with a computer program and thus obtain 
the partition function of the model exactly. On a con- 
ventional personal computer with a rather moderate pro- 
cessing time it is easy to obtain results up to Z = 13, and 
simulations |^ for larger polyions show that the qualita- 
tive behavior of the model does not change much beyond 
this value, so we will restrict ourselves here to Z < 13. 



III. RESULTS FOR THE THERMODYNAMIC 
BEHAVIOR OF THE MODEL 
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FIG. 2. Free energy ip as a, function of the angle 6 for Z = 9, 
n — 4, a = 0.8, C — 2. The curves shown are for: x = 2 (full 
line), X = 1.3676010335 (dotted line), and x = 1.3 (dashed 
line). Notice that the free energies for d = and — 7^/2 are 
equal in the second case. 

We will now focus on the regions of the parameter 
space where the parallel and perpendicular configurations 
are stable. We will start discussing the behavior of the 
model at vanishing values of x. In this limit, the parti- 
tion sum is dominated by the contributions coming from 
interactions between charges separated by a distance x. 
For the case 9 — we have Z pairs in this situation. We 
may then rewrite the partition function as 




The thermodynamic behavior of the model is deter- 
mined by the free energy cj), and we will start considering 
the dependence of (j) on 0. For convenience, we define an 
adimensional free energy 



Db 



\nQ. 



(13) 



The free energy is a function of the parameters Z, n, ^, 
a, X, and 9. For reasons which will become clear below, 
we replace the parameter a by 



2na 



1, 



(14) 



so that we will consider the free energy ip{Z, n, ^, a, x, 9). 
For all cases we noticed that the global minimum of the 
free energy is located either at the parallel {9 = 0) or at 
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V = minwi miny^ Vikk, 

i i ^ — ^ 

A;=l 



(15) 



where zq = Zkk — exp(— ^/x) vanishes as x ^ 0. Now let 
us call 



(16) 



and let us suppose that the first A^i of the Nc condensed 
counterion configurations correspond to this value of Vi. 
Thus 



Q|l =Zo^I^ll(l-KP), 



(17) 



where 
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Ni Z-1 



and 



P 



1 



=Afi + l 



n 



Z-l 



•Im 



n 



'Im 



vanishes as x — *■ 0, since Vi > v. The A^i configurations 
considered here are the ones that maximize the number 
of (— 1, — 1 + a) pairs separated by a distance x, and thus 
V = Z — 2na = —Za. For the paraUel configuration, we 
thus have that at small values of x, the partition function 
is asymptotically given by 



Q 



A simple combinatorial calculation leads to 

Z\ 



Ni = 



(n!)2(Z-2n)! 



(18) 



(19) 



The same line of reasoning may be applied to the per- 
pendicular case, where only the pair of central charges 
([Z+l]/2,l) and ([Z + l]/2,2) is separated by a distance 
X. One thus has in this case 




FIG. 3. Curves of equal free energies ((/jy = ip^) for the 
model with ^ = 9, n = 2, for ^ = 1 (full line), ^ = 2 (dashed 
line) and ^ — > oo (dotted line). The parallel configuration is 
stable below the curves. Notice, in the inset, the sign of ai 
varies with ^. in these cases. 



' 



(20) 



with 



^2 Z-l 

^-=Enj/r n 

i=l J = l (i,m)7^([Z+l]/2,[Z+l]/2) 



(21) 



where the first N2 configurations are now supposed to be 
the ones with a central (—1,-1 + a) pair, and 



No = 



{Z ~n)\{Z-n - \)\n\[n~ 1)!' 



(22) 



^From equations |ij and Q\\ — at small distances 
leads to 



a = oq + a\x + ©(x^), 



with 



Z-2n 
Z{2n- 1)' 
2n 

Zi[2n- 1) 



(lnW^2 -InVFi 



(23) 

(24) 
(25) 



It is thus apparent that the sign of the inclination of the 
curve = </3_L in the (a, a;) plane at a; = is determined 
by the sign of Wi — Wi, since ln(a;) is a monotonically 
increasing function. 



In Fig. [s] the curve of equal free energies in the plane 
(a, x) is shown for the case Z — 9, n — 2. For all tem- 
peratures, the curves meet at (oq = 5/27, a; = 0), as 
expected. For the case of the ground state (vanishing 
temperature), the curve displays discontinuous deriva- 
tives at points where the configurations of the ground 
state change. Another point which is worth observing 
is that for ^ > 1.489528... or ^ < 0.288876... the curve 
starts with negative inclination at (ao,0). Thus, for a 
value of a somewhat smaller than ao the perpendicular 
configuration is stable at large and small distances, while 
parallel rods are stable at intermediate distances. Thus, 
in this case we find that the perpendicular phase is reen- 
trant as the distance x is lowered. Figure || shows the 
initial inclination ai = (S^-) as a function of f for 
some examples. As is apparent in expression E5L these 
inclinations diverge as ^ ^ 0, their sign in this limit be- 
ing determined by W2 — W1. As a general rule, one notices 
that reentrant behavior is found for relatively small val- 
ues of n in the examples that we have studied (n < 3 for 
Z = 7,9, 11 and n < 4 for Z = 13). 
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Fig. H the curves /y = 0, /_i = 0, and = (p±_ are de- 
picted for a particular case. One notices that for a > oq, 
the force between the rods is attractive for sufficiently 
small distances x and the rods are in the parallel config- 
uration. For ax < a < ao the force is still attractive at 
small X, but the rods are in the perpendicular configura- 
tion. 



FIG. 4. Inclination of the curves of equal free energies, 
da/dx = ai, at a:: = 0, as a function of ^ for some cases 
{Z,n) from top to bottom: (13,4), (9,3), (11,3), (13,3), 
(9,2), (11,2). 

Finally, we discuss the behavior of the force defined in 
Equation 1^. As a general rule, for sufficiently large values 
of a, the force is repulsive at large distances and becomes 
attractive as the distance is lowered. It is useful to find 
the value of the distance x for which the forces /|| and /j_ 
vanish. The behavior of these curves at small distances 
may be found from the asymptotic behavior of the free 
energies described in equations |l^ and We find 



h 



Za 



(26) 



where 



1 dW\\ 



= h\\x + 0{x^). 



Thus, for small values of x, the curve /|| = reads 



(27) 



so that attractive forces in the parallel configuration are 
possible only if a > 0, since > 0. For the perpendicu- 
lar configuration, similar considerations lead to 



/j_ « h^x 



1 



(28) 



and therefore the curve = for x ^ 1 is given by 



a - 



1/3 



(29) 



For small separations between the rods, the curve /j^ = 0, 
in the (a,x) plane, tends to a = aj^ = —1 + 2n/Z. In 
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FIG. 5. Curves of equal free energies (full line), /y = 
(dotted line), and /x = (dashed line), for the model with 
Z = 9, n = 4, and ^ = 1. Note that bellow the solid curve 
the parallel configuration has the lowest free energy, while 
above the perpendicular configuration is energetically favored. 
The parallel configuration has an attractive force between the 
polyions bellow the dotted line. Similarly, the perpendicular 
configuration is attractive bellow the dashed line. 




FIG. 6. Force as a function of the distance between rods 
for Z = 9, n = 4, a = 0.2 and ^ = 1. The full line corresponds 
to perpendicular rods and the dashed line to parallel rods. 
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FIG. 7. Force as a function of me temperature 1/5 for 
Z = 9, n ^ A, a = 0.6, and x = 0.5. 

The behavior of the force as a function of the distance 
X is shown in Fig. ^ for a = 0.2. At a; = 1 the force is re- 
pulsive and the rods are perpendicular. As the distance is 
lowered, the force becomes attractive, with the rods still 
in the perpendicular configuration. At x = 0.246924... 
the rods change to the parallel configuration and the at- 
tractive force becomes much larger. At a small distance, 
the force displays an asymptotic behavior of the form 
— Za/x'^, according to Eq. At a = 0.75 the force 
is repulsive at large distances, with perpendicular rods. 
As the distance is lowered, the force changes discontinu- 
ously to attractive as the rods become parallel. Finally, 
an example of the behavior of the force as a function 
of the temperature 1/^ may be seen in Fig. |^. For the 
values of the parameters used to obtain these data, the 
parallel configuration is stable. In general, the modu- 
lus of the force increases as the temperature is lowered, 
since the charge correlations grow in this case. However, 
at relatively low temperatures and short separations be- 
tween the polyions, this rule may not apply, as is the 
case in the example shown. The reason for this is that 
the ground state configuration, which corresponds to the 
lowest electrostatic energy, is not, in general, the config- 
uration which maximizes the attractive force. For low 
temperature and short distances there are configurations 
which have forces more attractive than what is found in 
the ground state. The total force, being a weighted mean 
of forces associated with all configurations can, therefore, 
become more attractive than the force at T=0. 



divalent counterions. In the ground state the counterions 
are distributed periodically along the polyions. At finite 
temperature, this periodicity will be destroyed, the cor- 
relations between the condensed counterions, however, 
will persist. Thus, for sufficiently large electrostatic cou- 
pling, i.e. Manning parameter, the electrostatic potential 
at position (p, s) from a polyion can be approximated as, 



(-1)" 



^oo Vp^ + is - 



(30) 



Appealing to Poisson sum rule the asymptotic large dis- 
tance behavior or this sum can be evaluated yielding 



s) 



Aq 



: exp 



The energy of interaction between two lines of charge 
with staggered charges is then 



E 



and the force, F 




F 



Aq 



b 



(32) 



(33) 



We see that the correlation induced attraction decays ex- 
ponentially, with the characteristic range of Z = h/i:. 

We have studied the electrostatic interaction between 
two charged rods with a layer of condensed counterions. 
This is the simplest model of interaction for like-charged 
polyions in a polyelectrolyte solutions. It is shown that in 
spite of the equal net charge on the two macromolecules 
the correlations between the condensed counterions can 
produce an effective attraction. We show, however, that 
this attraction is of extremely short range, comparable 
to the monomer separation along the macromolecules. 
Furthermore, it is found that at large distances the 
monopolar repulsion between the charge densities forces 
the polyions into a perpendicular configuration. At short 
distances the correlations between the condensed counte- 
rions can become sufficient to produce a macromolecular 
alignment. The energy barrier associated with the transi- 
tion from the perpendicular to the parallel configuration 
might be relevant for the kinematics of bundle formation 
in solutions of stiff polyelectrolytes ||lj] . 



IV. DISCUSSION AND CONCLUSIONS 

To get a better understanding of the range of electro- 
static correlational forces involved, we study a simple ex- 
ample of two parallel long polyions with Z/2 condensed 
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